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Abstract: Heat transport in highly porous fiber networks is analyzed via two-point correlation functions. Fibers are 
assumed to be long and thin to for allow a large number of crossing points per fiber. The network is characterized 
by three parameters: the fiber aspect ratio, the porosity and the anisotropy of the structure. We show that the 
effective thermal conductivity of the system can be estimated from knowledge of the porosity and the correlation 
lengths of the correlation functions obtained from a fiber structure image. As an application, the effects of the fiber 
aspect ratio and the network anisotropy on the thermal conductivity is studied. 
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trinsic TC of the constituents and also on the specific 
topology of the fiber network and the matrix, which 
is a phenolic resin in the case of PICA. In fact when 
the TC of the matrix is negligible, the only conduc- 
tive heat transport mechanism in HPFM is the direct 
conduction through the fiber contacts. 

Several approaches and numerical methods have 
been reported in the literature for TC modeling of fi- 
brous material [3-8], White et. al., for example, inves- 
tigated numerically the influence of anisotropy on the 
electrical conductivity of networks of finite rods near 
the percolation threshold by using a resistor network 
(RM) approach [3]. In this approach, the fiber system 
is modeled as an equivalent network of resistors. The 
effective electrical conductivity was obtained by solv- 
ing a system of linear equations obeying Kirchoff’s 
laws at each node. It was shown that a discretized 
model can be used to calculate the ETC of fiber com- 
posites built from a low-conductivity matrix embed- 
ded with high-thermal conductivity fibers. Numerical 
results agreed with analytical models [6]. Qualitative 
estimation of conductivity has also been reported in 
Ref. [4]. It was shown there that the conductivity can 
be expressed in terms of the porosity, the fiber thick- 
ness and the mean distance between fiber crossings. 

Recently we showed that the ETC of HPFM can 
be estimated utilizing the two-point correlation func- 
tion (CF)[8], We demonstrated that the effective TC 
can be written in temis of the total porosity and the 
correlation lengths of the CF along each (x, y and z) 
direction. This approach is particularly useful for sys- 


1 Introduction 

Highly porous fibrous materials (HPFM) have many 
applications worldwide. Such materials range from 
paper products, flexible thin-film transistors, and 
physiological systems to high- temperature materials. 
In particular, they have extensive use in the aerospace 
industry as thermal protection systems. PICA (pheno- 
lic impregnated carbon ablator), a material developed 
by NASA, is an example of currently used HPFM- 
based insulating systems [1,2]. PICA is made from a 
carbon fiber insulation (Fiber Materials, Inc. under the 
trade name Fiberfom®) impregnated with a phenolic 
resin (for example, Borden Chemical SC 1008®). The 
typical diameter of the fibers in PICA is 14-16 //m 
and their length exceeds ~ 1 600 //m. The resin creates 
a highly porous themioset structure after polymeriza- 
tion with a low bulk density ranging from 0.22 to 0.27 
g/cm 3 . Characterization of the effective properties of 
this novel material is a challenging task. 

Themial conductivity (TC) is one of the most im- 
portant properties of HPFM. Since these materials are 
heterogeneous composites, the effective themial con- 
ductivity (ETC) is the relevant quantity. Accurate 
estimation of the ETC is important to provide input 
for themial response modeling under operating condi- 
tions and also to provide guidence for the selection of 
the optimum materials to meet specific design require- 
ments. In addition to low ETC, relevant properties 
of PICA include reasonable strength, low weight, and 
excellent ablative properties, among others. Themial 
conductivity in fibrous materials depends on the in- 



terns where the RM approach is computationally in- 
tractable. While that previous work presented a rigor- 
ous CF analysis of fiber networks and its implications 
for ETC determination, practical examples of the CF 
approach were not included. In the present work, we 
provide such examples and present results estimating 
the ETC of 2D models of fiber networks. We start 
in Section 2 by reviewing the main points of the pro- 
posed approach. In Section 3, numerical computa- 
tions of TC are performed and their dependence on 
the fiber aspect ratio and the anisotropy of the fiber 
ensemble is studied. 


2 Theory 

2.1 Thermal conductivity of fiber networks 

We consider a dilute random array of cylindrically 
shaped fibers in a matrix. We assume the matrix does 
not contribute to the thermal transport. Note we only 
consider conductive thermal transport, therefore, heat 
propagation through the fiber network is the only ther- 
mal transport mechanism in the structure. If n is the 
fiber concentration, i.e. the number of fiber in the unit 
cube, the volume fraction of the fiber phase can be 
written as Vf — imd?l/A (Vf — ndl in the 2D case), 
where l and d denote the length and the diameter of 
the fibers, respectively. The porosity of the structure 
is p — (1 — Vf ). We assume the fibers have large as- 
pect ratios, (l/d V 1), so that each of them has many 
crossing points along their length and thus the mean 
distance between the neighboring connecting points 
A (i.e. connectivity) is small compared to the fiber 
length (A <C 0- This also means that the fiber concen- 
tration is large enough, so that we are well above the 
percolation threshold. 

We also assume cylindrical symmetry for the fiber 
network along the 2 -axis and consider an external 
thermal gradient applied along this axis, i.e., the heat 
flow takes place along the 2 -direction. The heat flux, 
Iij, through a fiber connecting two points r, and r ;l 
at two close planes z = 0 and z — h, respectively is 
given by 


h 3 = r ^-[T{r l )-T{r ] )}. 


Here T(r,;) denotes the temperature at ry, t;j is 
the fiber thermal conductivity, l r] = \r l: — r, = 
h/\ cos Oij |, and Oij is the angle between the fiber and 
the 2 -axis. The total flux through the plane z = 0 
is obtained by adding all fluxes corresponding to the 
fibers crossing the plane 
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This allows us to determine the effective thermal 
conductivity as: 


ft'eff = 
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In the above expression, N — AS/tiD\ is the total 
number of fibers crossing the plane at 2 = 0, S' is the 
side area of the specimen in this plane, and D± is an 
average distance between the fibers crossing the plane 
2 = 0. Supposing that (cos 6V [T(r,;) — T(r,-)1) = 
(cos 0ij){[T(ri) - T(rj)\), one finds 
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or 


«eff = 
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where D\\ is an average distance between the fibers 
along the 2 -direction, so that Vf = d 2 (2Dj_ + 
D\\)/{DlD{). 

Following a similar analysis we can consider the 
TC across the 2 -axis. The results for the components 
of « e ff along and across the the 2 -direction are respec- 
tively 
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Equations (1) and (2) express the fact that the ef- 
fective TC of the fiber network can be estimated by 
knowing only: (i) the intrinsic TC of the fiber, (ii) 
the components of the mean distance among the fibers 
along and across the heat flux direction and (iii) the 
porosity of the structure. Such expressions can also 
be used for 2D fiber networks by omitting the factor 2 
in the denominator. In the next section, we show how 
the D parameter can be estimated from the CF of the 
fiber structure. 


2.2 Two-point correlation function 

Valuable information on the morphology of heteroge- 
neous materials is contained in the the n-point corre- 
lation functions (n — 1, 2, 3, • • • ). These microstruc- 
tural descriptors were originally introduced by Brown 
[9] in the context of the effective properties of com- 
posites and subsequently have been deeply investi- 
gated by Torquato [10], 

We introduce the two-point CF of the fiber struc- 
ture as 


w(r /x = (b?(r) - (t?(r))][?7(rQ - (??(/))]) 
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Figure 1: Correlation function W{x, 0) of a 2D fiber 
network (b), measured along the horizontal axis of the 
image, x. All fibers have the same length, l — 20 
units and width d — 0.4 units. Fiber centers and ori- 
entations were allow to have random values. Typical 
correlation lengths of W{x, 0) are highlighted in (a). 

where (• • • } denote an ensemble average and the char- 
acteristic function i/fr) given by 

. . f 1. if r is inside a fiber 

nlr) — < 

[ —1, otherwise, 

represents the fluctuation of the local TC with respect 
to its mean value. For statistically homogeneous me- 
dia, equation (3) depends only on the distance be- 
tween two arbitrary points r and r ' , [W{r, r') — 
W{\r — r'|)l; it is equal to unity at the coordinates 
origin, W{0) — 1, and vanishes at the infinity (see 
Fig. 1). 

The CF as defined by equation (3) can be directly 
computed from the digital image of the fiber structure. 
There currently exists a variety of nondestructive tech- 
niques to obtain two- and three-dimensional images of 
materials, including optical microscopy, synchrotron- 
based tomography and magnetic resonance imaging 
[10]. Through a thresholding procedure, the gray- 
scale image of a fibrous material can be reduced to 
a binary image in which gray values lighter than a 
chosen threshold are set to white and the others set 
to black [11]. The resulting pixel data can then be 
described by the stochastic function, //(r), and equa- 
tion (3) can be directly applied. 

The profile of the CF curve provides useful in- 
formation about the fiber structure. For example, the 
fiber length (i.e., largest fiber dimension) can be esti- 


mated as the distance at which the CF drops to zero. 
For the fiber structure illustrated in Fig. 1, this value 
is found to be If — 50 d, which is in perfect agreement 
with the fiber length used in the simulation where we 
had set If — 20 units and d — 0.4 units. The first cor- 
relation length of the CF curve, 1 C1 , on the other hand, 
is related to the smallest fiber dimension: the width 
of the fibers. One interesting question is: to which 
parameter of the fiber structure is associated the sec- 
ond correlation length (second inflection point), l C2 , 
appearing in the inset of Figure 1? 

As demonstrated in [8], this parameter provides 
information about the fiber network as a whole. For 
2D fiber structures, where all fibers intersect if they 
are long enough and close, the second correlation 
length coincides with the fiber connectivity (average 
distance between interfiber crossings), l C2 ~ A, which 
determines the conduction pathways in the structure. 
For 3D fibers structures on the other hand, the second 
correlation length is not A but the mean distance be- 
tween the fibers, l C2 ~ D. Therefore A or D can be 
estimated from equation (3) and inserted into equa- 
tions (1) and (2) to estimate the effective TC. 


3 Numerical Results and Discussion 

In this section, we present results that used the ap- 
proach presented above to estimate the effective TC 
of 2D fiber structures as a function of the fiber param- 
eters. The structures simulated were constructed of 
conducting sticks (i.e., straight fibers) of equal length 
l and width d (l d) and were randomly positioned 
in a unit-size square. Within this square we generated 
random ( Xi , yi),i — 1,2 ■ ■ ■ N, pairs of coordinates 
by using a pseudorandom number generator. We used 
the same “seed” to allow the same data to be gener- 
ated for repeated evaluations. We then attached the 
center of one stick to the site ( Xi,yi ) by choosing a 
random orientation 0 t , with respect to the y-axis com- 
prised in a given interval, — 0l < 0, < 0 L . For these 
values, a stick was plotted by drawing a rectangle be- 
tween the points [ Xi — {1/2) sm0i,yi — {1/2) cos 0/\ 
and [xi + {1/2) sindi, yi + {1/2) cos 0/. By changing 
the limiting values, ±9l, different degrees of align- 
ment for the fibers can be obtained and thus different 
anisotropies for the fiber ensemble produced. The de- 
gree of anisotropy is quantified by an orientation pa- 
rameter, S — 2 (cos 2 0/) — 1, which varies from, 0, for 
isotropically oriented sticks, to 1, for perfectly aligned 
sticks. The fiber network is thus fully characterized by 
three parameters, namely the porosity, p, the fiber as- 
pect ratio, (l/d) and the anisotropy, S. Different val- 
ues of these parameters lead to different fiber struc- 
tures. Figures 2(a) and 2(c) show two 1000-sticks 



Figure 2: (a) Two dimensional fiber network with porosity p — 0.81 and orientation parameter S — 0.01, composed 
by 1000 fibers with aspect ratio l/d = 100. (b) Correlation functions W(x, 0) and W(0,y) of (a) taken along the 
horizontal and vertical axis, respectively, (c) Same as (a) but for S — 0.40. (d) Same as (b) but for (c). 


samples with fixed (l/d) — 100 and p — 0.81 for 
S = 0.01 and S = 0.40, respectively. 

We started by analyzing the effect of the fiber as- 
pect ratio on the effective TC. To do that we fixed the 
width of the fibers to a length of 0.0015 unit (i.e., 
d — 1.5% of the square side length) and we chose 
an anisotropic fiber configuration characterized by a 
S — 0.4 value. A total of 50 structures with different 
aspect ratios in the interval, 80 < (l/d) < 2000, were 
generated to form high-resolution images and used as 
input for TC estimation. A full data set of 1400 x 1400 
pixels for the images was enough in providing good 
statistic for the CF computations. In our analysis we 


assumed that the heat flow in the structures take place 
along the x-direction (i.e., horizontal axis in Fig. 2). 

Following the procedure outlined in Section (2.2), 
the CF of each fiber network was calculated along and 
across the heat flow direction and the second correla- 
tion length, l C2 , was estimated from them. The fiber 
density Vf was evaluated by counting the number of 
pixels belonging to the fiber phase (“black pixels”). 
After inserting values of fiber density and correlation 
length into the 2D version of equation (2), the ratio 
of the effective TC to the intrinsic fiber conductiv- 
ity, K e f i/nf, was calculated for each simulated struc- 
ture. Figure 3(a) presents the results as a function of 
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Figure 3: Effective thermal conductivity along (/c r ff ) and across (k% s ) the normal to the substrate (coincides with 
direction of heat flow) as a function of the fiber aspect ratio (a) and the orientation parameter S (b) (see text for 
details). 


l/d. As it can be seen, the effective TC along and 
across the heat flow direction, («* ff , n y e ^), increases 
with increasing fiber aspect ratio. The component 
of TC along the heat flow, for example, starts with 
a low value of i; x lf = O.OG/Vf when the percolated 
network is formed by “short" fibers (l/d ~ 80), and 
then increases slowly until it reaches a constant value 
of about 0.28h'f . A similar behavior is also observed 
for K^ ff , but in this case, TC values saturate at large 
l/d to 0.48/, y. This is an expected result, since as 
the fibers become longer, the number of fiber con- 
tacts increases and therefore the connectivity value of 
the fiber network decreases, giving an increase in TC. 
When the fibers are long enough, further elongation 
of the fibers means an increase in the number of dead 
ends. However, this does not change the connectivity 
of the structure and the effective TC remains constant. 
This is precisely the behavior observed for both n x K 
and K y s at large aspect ratio values, l/d> 1300. 

The gap observed in Figure 3(a) between the /v* f 
and n y R curves is a consequence of the anisotropy of 
the structure, as we used a fixed value of S — 0.4 for 
all the simulations. Thus, for each sample, the fibers 
were generated with a significant degree of align- 
ment along the y-axis [see Figure 2(c)]. As a con- 
sequence, the l C2 of the CFs along the //-axis was 
found to be larger than those measured along the heat 
flow direction [see Fig. 2(d)]: therefore, according to 
equation(2), /^ ff > «® ff ,V(Z/d). 

Next we studied the influence of anisotropy on the 
effective TC. In order to do that, we simulated fiber 
samples with the same aspect ratio (l/d — 100) and 
porosity (p = 0.81), but with different orientation pa- 


rameter values. The estimation of the effective TC 
was carried out analogously to our previous study, that 
is, through the CF analysis of the respective fiber im- 
ages. The results are shown in Figure 3(b). 

From the figure it can be seen that there are two 
characteristic trends in the effective TC as a function 
of the network orientation. First, for low values of S, 
the effective TC along and across the heat flow direc- 
tion are roughly the same, = K y s ~ 0.099/if. 
This is not a surprising result, since at those values 
of S the fiber network is strongly isotropic and there- 
fore the correlation lengths of the CF along any di- 
rection coincide [see Fig. 2(a)], Second, as the axial 
fiber alignment (S > 0.1) increases, H x tf and n y s de- 
part from 0.099/if in a completely different way: /c“ ff 
drops to 0.002/tf while n y s increases until reaching 
a values of 0.2/if at S — 1. The explanation of the 
above relies on the fact that as the fibers become more 
favorably oriented along the j/-axis, the number fiber 
contacts along this axis decreases while the contacts 
along x increase. This change in the network architec- 
ture causes an increase in connectivity with S along y 
and thus a decreasing (increasing) of the effective TC 
along (across) the heat flow direction. The same de- 
pendence of Kgg for 3D structures was found by White 
and coworkers [3], 

4 Conclusion 

We have shown that the effective TC of percolated 
fiber networks can be estimated by utilizing informa- 
tion contained in the CF of the structure. The only pa- 
rameters required are the porosity of the structure (or 


the fiber density) and the correlation length of the CF [11] J. G. Berryman, Measurements of spatial cor- 
along some characteristic directions. The estimations relation functions using image processing tech- 

of TC shown in this work for 2D networks demon- niques, J. Appl. Phys. 57, 1985, pp. 2374-2384. 

strate the applicability of the CF approach. 
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